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Abstract. We suggest an equal-time n-point correlation function for systems in 
the directed percolation universality class which is well defined in all phases and 
independent of initial conditions. It is defined as the probability that all points are 
connected with a common ancestor in the past by directed paths. 



1. Introduction 

One of the most studied classes of continuous phase transitions far from equilibrium is 
directed percolation (DP) pQ. This transition occurs e.g. in stochastic lattice models 
for epidemic spreading, where infected sites either recover or infect a randomly selected 
nearest neighbor. Depending on the rates for infection and recovery, an initial epidemics 
may spread or disappear. Once the infection becomes extinct, the system enters a 
so-called absorbing state from where it cannot escape. Both regimes, survival and 
extinction, are separated by a well-defined phase transition. Directed percolation is the 
simplest universality class of phase transitions into absorbing states [2H7] and plays a 
paradigmatic role similar to the Ising model in equilibrium statistical mechanics. 

Critical phenomena can be studied by analyzing correlation functions. In DP an 
important example is the so-called pair- connectedness function c(fi, t\\ r*2, £2) which is 
defined as the probability to find a directed connected path from a lattice site located at 
position fi and time t\ to another site at position r*2 and time ti > t\ (see e.g. [7]). The 
pair connectedness function is well defined below, at, and above the critical point and 
plays the role of a propagator in the field-theoretic formulation of DP. However, if one 
is interested in spatial correlations at a given instance of time the pair connectedness 
function does not provide any information because it vanishes for t\ = £2- This raises 
the question how a meaningful equal-time correlation functions can be defined in models 
with absorbing states. 
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Figure 1. Equal-time two-point correlation functions. Left: Conventional density- 
density correlation function in 1+1 dimensions, defined as the probability that there arc 
two directed paths connecting the initial state with the two points. Right: Correlation 
function investigated in the present work, defined as the probability that there exists 
a common ancestor in the past which is connected with both points by directed paths. 



The simplest equal-time correlation function, which has been studied in various 
contexts, is the density- density correlation function (s(fi, t)s(r 2 , t)), defined as the 
probability to find two sites at positions fi and r 2 simultaneously in the active (infected) 
state at time t (see Fig.[T]). This two-point function seems to be a natural candidate since 
it reminds us of a spin-spin correlation function in equilibrium statistical mechanics. 
Unfortunately this correlator depends on the history of the process. For example, if 
the system has not yet reached a stationary state, its value will depend on t and the 
specific choice of the initial state. As a way out, one could measure the correlations in 
the stationary state, but then one is restricted to the active phase and the result will 
depend on the distance from criticality A = p — p c . To study the critical properties of 
the correlator, one has to keep the system in a slightly supercritical stationary state 
< A <C 1 and to restrict the analysis to distances below the correlation length 
r <C r c ~ A~ u± . This method is of course numerically inefficient and is further 
complicated by crossover effects. 

To overcome this problem, Dickman et. al. [8] devised a method that allows one to 
study equal-time correlations at the critical point. Starting point is the observation that 
the density of active sites in a finite system of lateral size L at criticality first decays as 
usual, then stays for a while in a quasi- stationary state with a density p ~ I//"-*-, until it 
suddenly enters the absorbing state. The key idea is to prevent the system from entering 
the absorbing state. This is done by keeping a list of typical configurations sampled in 
the quasi-stationary state and to reimpose one of them whenever the system would 
have entered the absorbing state, thereby artificially keeping it in its lowest excited 
state. Although this method became one of the standard tools in the study of absorbing 
phase transition, it is not entirely clear how the artificially stabilized state has to be 
interpreted. Moreover, the method involves finite-size effects which have to be handled 
with care. 

In this paper we propose a different type of equal-time correlation function which 
circumvents these difficulties. This correlator is defined as the probability that two 
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points located at r\ and r*2 at equal time have a common ancestor somewhere in the 
past, meaning that this ancestor is connected to both points by directed paths (see right 
panel of Fig. [1]). Note that it only matters whether this ancestor exists, but not where 
it is located in space and time. Obviously, this definition can be generalized easily to 
n-point functions. 

Compared to the conventional density- density function, the correlator proposed 
here has several advantages. As it is defined in terms of paths, it is by definition 
history- independent and well-defined in the limit L — > oo. Moreover, it can be analyzed 
everywhere in the phase diagram so that the trick by Dickman et. al. is no longer 
needed. 

The main motivation for the present work came actually from a different direction. 
In equilibrium statistical mechanics the critical systems are usually not only scale- 
invariant but also conformally invariant. Especially in two dimensions this symmetry 
restricts the possible universality classes strongly, leading to a finite classification scheme 
of such transitions. There have been various attempts to generalize this concept to non- 
equilibrium systems by introducing an anisotropic version of conformal invariance (see 
e.g. [9] and references therein). This concept was successfully applied to a large number 
of models but seemed to fail in the case of DP [TU|[TT]. The ongoing debate triggered the 
question whether time-independent correlations of a two-dimensional DP process would 
exhibit not only scale-invariant but also conformal invariant features. This motivated 
us to search for a well-defined equal-time n-point function in DP, where the predictions 
of conformal invariance can be tested. At first glance it turned out that the correlation 
function defined above behaves like a conformally invariant one. However, extensive 
simulations reveal small violations. 

2. Numerical implementation 

The equal-time n-point function $(fi, . . . , f n ) introduced in this paper is defined as the 
probability that there exists a lattice site in the past, a so-called common ancestor, 
which is connected to all points by directed paths. 

Let us first discuss the question how this correlator can be estimated efficiently in 
a numerical simulation. To this end one has to choose a realization of DP on a lattice. 
For reasons to be explained below, we chose directed bond percolation, where the bonds 
of a diagonal lattice are open with probability p and closed otherwise (see e.g. (2HZ])- 
Activity percolates through open bonds along a preferred direction denoted as time. 

To estimate the correlator defined above, we could generate a spatio-temporal 
configuration of open and closed bonds, select n points with given distances and trace 
along the open paths backwards in time to see whether all points have a common 
ancestor. This is shown in the left panel of Fig. [2] in the example of a three-point 
function in 1+1 dimensions. However, since the ancestor can be located anywhere, no 
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Figure 2. Example of a three-point function in 1+1 dimensions. Left: Configuration 
of open and closed bonds, where the three points marked by red dots have a common 
ancestor. Middle: The same configuration turned upside down may be interpreted as 
a time-reversed DP process starting with three active sites. The point where the three 
generated clusters merge for the first time corresponds to the common ancestor in the 
original setting. Right: The merging point is determined most easily by percolating 
bit patters (see text). 

matter how far it goes back time, any spatial or temporal restriction of the generated 
configuration of open and closed bonds will inevitably lead to finite-size effects. 

A more efficient method exploits the so-called rapidity reversal symmetry of DP 
which states that directed paths forward and backward in time have the same statistical 
weight. This allows us to reverse the direction of time, use the n points as initial seeds 
of activity, and search for a common descendant in the future. Note that the rapidity 
reversal symmetry is exactly fulfilled in the case of directed bond percolation, which 
explains why we chose this particular realization of DP. 

To find out whether n points have a common descendant, imagine that each seed 
carries a color, as shown in the central panel of Fig. [2j By definition the common 
descendant is a lattice site where all these colors mix for the first time. In a numerical 
simulation this intuitive color scheme can be implemented by associating a string of n 
bits with each lattice site and assigning one of these bits to each of the seeds (see right 
panel of Fig. [2]). These bit patterns then percolate forward in time along the open bonds. 
If two different patterns arrive at the same site, they coagulate by an OR operation. 
The first site carrying the bit pattern 11 ... 1 is a common descendant, corresponding to 
the common ancestor in the original situation. In this case the subroutine would return 
a '1', indicating that a common ancestor exists. However, if one of the bits dies out 
before a common descendant has been identified, the simulation stops and returns a '0'. 
The actual value of the correlator is then estimated by averaging over many independent 
runs. 

Technically, instead of using static arrays, it is advantageous to keep the coordinates 
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d 


P 


v± 




Z = V\\/u±_ 


X = P/v± 


6 = p/v\\ 


1 


0.276486(8) 


1.096854(4) 


1.733847(6) 


1.580745(10) 


0.252072(12) 


0.159464(6) 


2 


0.5834(30) 


0.7333(75) 


1.2950(60) 


1.7660(16) 


0.7955(30) 


0.4505(10) 



Table 1. Numerical estimates of the DP critical exponents in d+l dimensions, see [7] 
and references therein. 



of sites with a non-vanishing bit pattern in a dynamically generated list without 
restricting the range of the coordinates. This technique accelarates the simultion 
significantly and eliminates possible finite-size effects. Moreover, it is useful to perform 
the simulation without a cutoff time until a '0' or a '1' is returned, which eliminates 
possible finite-time effects as well. 

3. Scaling properties 

The phenomenological scaling theory of absorbing phase transitions involves four 
independent critical exponents /3, 0', u±, v\\ for the order parameter and its reponse 
field as well as the spatial and temporal correlation lengths [7]. In DP, the rapidity 
reversal symmetry implies that (3 = (3', leaving three of the bulk exponents independent. 
Moreover, it is useful to define the derived exponents 

X = P/v ± , 8 = P/v\\ z = u\ { /u ± . (1) 

In d < 4 spatial dimensions the values of the critical exponents are not known exactly. 
Numerical estimates are listed in Table [3J 

Generally, an n-point correlation function involves n fields, each of them with a 
specific scaling dimension. In DP there are two primary fields, namely for creating 
and measuring activity. Because of the rapidity reversal symmetry they carry the same 
scaling dimension [x]~ x , where [x] denotes the unit of length. For example, the pair 
connectedness function c(r, t) mentioned in the introduction involves two fields and 
therefore carries the dimension [x] _2x (or equivalently [t]~ 2<5 , where [t] = [x] z denotes 
the unit of time). This implies that for r = the pair connectedness function has to 
decay as c(0, t) ~ t~ 2S . 

The correlation function $(r 1; . . . ,r n ) considered in the present work depends on 
n points and involves the common ancestor as an additional point. However, it is 
reasonable to assume that this point does not contribute to the scaling dimension of 
the correlator because its location in space and time is arbitrary, meaning that it is 
effectively integrated out. Therefore, the scaling dimension of this correlator is [x] _nx . 
For example, the two-point function is expected to decay as 

®{ri,r 2 ) ^ C\ri - r 2 \~ 2x , (2) 
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Figure 3. A: Numerically determined two-point function $(r) in 1+1 dimensions 
(black) compared with the expected power law r~ 2x (red). B: Same data multiplied 
by r~ 2x . C: Low distance correction (see Sect. [3]). 



where C is a non-universal amplitude. As shown in Fig. [3J, this power-law decay is 
accurately reproduced in numerical simulations with C = 1.145(25). 

The three-point function $(?+, r 2 , r 3 ) carries the dimension [x] _3x and depends on 
three distances ?"i2, ri 3 , r 2 3, where r^ = jr, — fj\. Therefore, it is expected to scale as 

$(n,r 2 ,r 3 ) ~ (r 12 r 13 r 2 3)~ x F (— , — J (3) 

with a universal scaling function F that depends only on scale-invariant ratios of these 
distances. 

As usual, scaling laws are only valid in the asymptotic limit of large distances, 
while for short distances one expects non-universal corrections to occur. For example, 
the two-point function shown in Fig. [3] clearly deviates from the predicted power law 
over the first two decades of time. In fact, the correlator is predicted to diverge if 
r = |r*i — r*2 1 goes to zero. On the lattice, however, $(0) = 1 is always finite. 

We propose that these short- distance corrections can be taken into account in a 
lowest-order approximation by shifting the distance parameter by a constant, i.e. Eq. (J2J) 
is replaced by 

$(r) ~ C(r + r y 2x . (4) 

where r is a constant. This constant should be chosen in such a way that the short- 
time deviations are minimized, as demonstrated in Fig. |3p, leading us to the estimates 
r ~ 1.05(10) in d = 1 and r ~ 0.31(1) in d = 2 spatial dimensions. Note that these 
values are expected to be non-universal, i.e. they may differ for different realizations 
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Figure 4. Three-point correlation function in d = 1 spatial dimension. The graphs 
show the value of the scaling function -^(77^ 77^) — &{ r ii r 2, ? , 3)/(?'i2'"i3?'23) x for 
ri3 = 256 fixed and ri2 varying between and 128. The red dots with error bars 
report two precision measurements which confirm the drift. 

of DP. Moreover, we would like to point out that the constant should not be derived 
from the condition $(0) = 1 since in this point higher-order corrections would lead to 
an incorrect estimate. 

4. Test of conformal invar iance 

In 1970, Polyakov argued that scale-invariant equilibrium systems should also be 
invariant under local scale transformations [12]. This conformal invariance greatly 
restricts the form of correlation functions. In particular Polyakov showed that the 
three-point function at criticality should take the form 

G IH {ri, f 2 , f 3 ) = const rf 2 3 " Al ~ A2 rf^ 1 "^ rf^ 2 ^ 3 , (5) 

where = |fj — r}| and Aj is the scaling dimension of the field measured at position r*j. 
This equation is expected to hold in any dimension. Later it was shown that conformal 
invariance in two dimensions, where the conformal group has infinitely many generators, 
is even more restrictive and leads to a complete classification scheme of equilibrium 
critical phenomena [T3j[T4]. 

In non-equilibrium statistical mechanics such a classification scheme is not yet 
known. Although there have been attempts to generalize conformal invariance to 
anisotropic and time-dependent systems [9], DP turned out to be incompatible with 
the predictions of such theories pUHH] . This led to the question of whether at least the 
time-independent properties of DP carry a signature of conformal invariance, motivating 
us to define a history-independent equal-time correlator. According to Polyakov, the 
three-point correlator defined above with scaling dimensions Ai = A2 = A3 = x) is 
conformally invariant if 
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Figure 5. Scaling function F(^, of the three-point function in 2+1-dimensions. 
The distances 7*12, 7*13, and vary between 5 and 90 lattice sites. The results are 
shown on the left side in a one-dimensional projection by plotting F as a function of 
TnTizTTi- The right panel shows the same data including short-distance corrections 
according to Eq. (jl]). The error bars refer to the statistical error only, they do not 
include possible systematic errors. 
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meaning that the scaling function F in Eq. (E]) is constant. 

To test this prediction, we first consider the one-dimensional case, where the three 
points ?r, 7^2,^3 lie on a straight line. In this case ^ + ^ = 1 so that the scaling 
function depends only on one argument, say The numerical results are shown in 
Fig. HI Plotted on a scale between and 1 (left panel), the scalling function seems to 
be almost constant. However, the magnification of the data shown in the right panel 
exhibits a clear systematic drift of roughly one percent. High-precision simulations at 
individual points confirm that this drift does not go away as one increases the distances 
even further, meaning that it cannot be attributed to lattice- or finite-size effects. 

To investigate the correlation function in 2+1 dimensions we performed simulations 
using a fixed lattice with 1024 x 1024 sites and periodic boundary conditions. In order 
to increase the numerical efficiency, we implemented a parallel version of the algorithm 
based on the programming language CUDA developed by NVidia, which allows one to 
access the multi-processor GPU of a graphics adapter. The point t\ was fixed during 
the simulation, while the locations of the second and the third point varied in the two- 
dimensional plane within the radius of 5 (border to the higher-order lattice effects ) to 
90 lattice sites from the first point. At these relatively small distances we expect finite 
size effects to be negligible. 

The results shown in Fig. [5] confirm the findings in one dimension: The scaling 
function is finite and has no signatures of zeros or singularities. It varies more strongly 
than in the previous case but still in a range of less than ±10%. 



Equal-time correlation function for directed percolation 



9 




X 



Figure 6. Left: Anisotropic limit of a three-point function in the time-reversed setting. 
In practically all configurations, where the three points have a common descendant, the 
colors of the nearby points mix already after very short time, thus effectively behaving 
like a single seed (right panel). 



5. Strongly Asymmetric limit 

In statistical physics most scaling function show a singular behavior, i.e. they vanish or 
diverge if one their parameters goes to zero or infinity. The scaling function F, however, 
remains finite even when two of the three points approach each other. As we will show 
in the following, the correlation function can be factorized in this limit. 

To understand this factorization intuitively let us consider two clusters generated 
two nearby points in the time- reversed picture (e.g. the red and the green one in Fig. |6]). 
Generating such clusters and selecting those which survive for a long time one observes 
that almost all of them consist mainly of a single color (red, green or its mixture), while 
the coexistence of different colors is exponentially suppressed for large t. Only clusters 
surviving with a mixture (symbolized by the orange line in Fig. E]) can contribute to the 
three-point function if they merge with the cluster generated by the distant point. 

From the perspective of the common descendant the points with small distance on 
the left side behave effectively like single seed (see right panel) with a certain survival 
probability. Knowing this survival probability and the two-point function, this allows 
one to predict the three-point function in the limit that two of the three points are very 
close. 

More specifically, let r = r 12 be the small distance and R ~ r 13 ps r 23 the mean of 
the large distances between the particles so that we can express the three-point function 
<&(R, r) as a function of these two distances. Moreover, let P\(t) be the usual survival 
probability of a cluster starting with a single seed and P^ir, t) be the survival probability 
that a cluster starting with two seeds of different color at distance r will survive in a 
mixed-colored state. Dimensional power counting implies the scaling laws 



with certain non-universal amplitudes A and B. The factorization conjecture is based 
on the assumption that for r <C R the three-point function $(i?, r) is proportional to 
the two-point function $(i?) because the pair of nearby particles can be replaced by an 





(7) 



Equal-time correlation function for directed percolation 



10 
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A 


B 


C 


F(0,1) = ^ 


1 


impossible 


0.968(5) 


0.9325(5) 


0.78(1) 


1.145(25) 


0.957(10) 


2 


0.334(2) 


0.294(2) 


0.9464(16) 


0.36(1) 


0.752(2) 


0.286(9) 



Table 2. Special values of the scaling function ([3]) for an equilateral triangle 
ri2 = ri3 = r23, three points on a line with equal distances r\i = r23 = and 
in the limit r 12 <C ?'i3 ~ ^23- 



effective single seed. However, the survival probability of a pair in a mixed state differs 
from the ordinary survival probability of a single seed, meaning that the proportionality 
factor is just P2(r,t)/Pi(t), hence 

P (r t) BC 
<$>(R,r) ~ 2 y/ $(ifl ~ —r x R- 2x for r<i?. (8) 

Since the time dependence cancels out, this result is compatible with the scaling form ([3]) 
with the special value of the scaling function 

F(0,1) = ^. (9) 

Determining the constants A,B,C numerically the values of this expression are 
compatible with the values obtained in direct numerical simulations (see Table [2]). This 
confirms that the scaling function does not diverge, even if two of the three points 
approach each other. Moreover, the numerical values support the observation that the 
scaling function is not constant but varies in a narrow range. 

In Table [2] we can see, that an equilateral triangle provides the higher value of F, 
than a straight line. Considering other constellations we saw that all other values of F 
vary between these two. We plotted F as a function of so-called normalized area of the 
triangle 

A N = A triang i e / (r 12 r23r 13 ) 2/3 (10) 

and could recognize a certain correlation between these variables (see Fig. [7]). Both of 
them do not exceed a certain range: F G [0.294(2); 0.334(2)] and A N G [0, a/3/4], and 
F grows with A^. Although it was not possible to determine an exact function due 
to a large error range, this observation can be helpful for further study of the scaling 
function. 



6. Conclusions 

In this work we have introduced a novel kind of equal-time correlation function for 
directed percolation which is defined as the probability that n points have a common 
ancestor in the past, no matter where this ancestor is located. This correlation function 
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Figure 7. Scaling function F with short-range corrections as a function of the 
normalized area Am (see Eq. (|TU)) ). Am = means that the three points lie on a 
straight line, and Am = 0.433 (maximum value) represents an equilateral triangle. 



is well-defined in the whole phase diagram and does not depend on the history of the 
process. Using the rapidity reversal symmetry, this correlation function can be related 
to the probability that the clusters generated by n seeds will merge at some point in 
the future, which allows one to compute the correlator efficiently 

Analyzing the scaling behavior of the three-point function we find that the universal 
scaling function is almost constant, which would be an indicator of conformal invariance. 
However, extensive simulations reveal a small variation of (in one dimension) less than 
one percent. 

The bottom line of this paper is a warning: Some properties of directed percolation 
behave almost as if the process was conformally invariant, but high-precision simulations 
reveal that this is not the case. This means that numerical results can easily generate 
misleading expectations and have to be interpreted with care. 
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